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LONG-TERM  GOAL 

The  long-term  goal  of  this  project  is  to  develop  computationally  efficient  and  optimal  techniques  for 
the  direct  assimilation  of  satellite  altimetry  data  into  the  global  and  mesoscale  ocean  circulation  models 
in  order  to  accurately  reconstruct  the  ocean  circulation  fields  thereby  improving  the  capability  of  ocean 
nowcast  and  forecast. 

OBJECTIVES 

The  primary  objective  of  our  work  is  to  interface  fast  implementations  of  the  Kalman-Bucy  filter  (KBf) 
that  we  have  developed  [1,  2,  3,  4,  5]  with  the  Naval  Research  Laboratory  (NRL)  Layered  Ocean 
Model  (NLOM)  and  then  carry  out  extensive  data  assimilation  studies.  We  are  working  in  close 
collaboration  with  the  group  of  Dr.  Harley  Hurlburt  at  the  NRL  Stennis  Space  Center. 

APPROACH 

We  assimilate  altimetry  data  into  NLOM  by  decoupling  the  ocean  circulation  fields  into  deterministic 
and  random  components.  We  apply  dynamic  linearization.  Dynamic  linearization  involves  two  stages 
at  each  time  step:  the  nonlinear  deterministic  stage,  and  the  linear  stochastic  stage.  With  the  nonlinear 
deterministic  stage,  the  nonlinear  stochastic  primitive  equations  is  integrated  forward  during  the 
current  time  step  from  assumed  deterministic  known  boundary  conditions,  known  initial  conditions, 
and  known  driving  forcings.  This  first  nonlinear  stage  uses  the  full  nonlinear  set  of  primitive  equations 
but  considers  only  the  nominal  operating  conditions:  boundary  conditions,  initial  conditions,  and 
forcings.  The  second  linear  stochastic  stage  of  dynamic  linearization  linearizes  the  nonlinear  primitive 
equations  about  the  deterministic  nonlinear  solution  obtained  from  the  first  stage.  The  nonlinear 
deterministic  equations  in  the  nonlinear  stage  of  dynamic  linearization  are  obtained  with  the  NLOM. 
The  linearization  of  the  primitive  equations  in  NLOM  gives  the  state  equation  for  the  KBf.  The 
measurement  equation  for  the  KBf  is  provided  by  the  satellite  altimetry  scan  pattern.  The  KBf  can  then 
use  the  linearized  equations  and  the  satellite  measurements  to  estimate  the  random  components  of  the 
circulation  fields.  Once  these  estimates  are  available,  they  are  used  to  correct  at  the  current  time  step 
the  nonlinear  deterministic  field  components.  These  corrected  nonlinear  deterministic  components 
provide  then  the  initial  and  boundary  conditions  for  the  nonlinear  stage  of  dynamic  linearization.  The 
process  repeats  at  every  time  step.  When  no  new  measurements  are  available,  the  dynamic  linearization 
keeps  marching  forward  through  the  nonlinear  primitive  equations. 
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Direct  implementation  of  the  KBf  leads  to  a  computationally  unfeasible  problem.  To  overcome  this  we 
exploit  the  block  structure  of  the  model.  In  particular,  the  system  matrices  associated  with  the 
linearized  set  of  primitive  equations  and  the  KBf  exhibit  a  banded-block  structure  that  enables  us  to 
reduce  significantly  the  computational  cost  of  the  KBf.  The  local  nature  resulting  from  the 
discretization  of  the  state  and  measurement  equations  of  the  KBf  is  naturally  reflected  in  the  block- 
banded  structure  of  the  system  matrices,  with  blocks  themselves  being  banded  and  sparse.  By  coupling 
this  block  structure  with  the  sparsity  of  the  satellite  measurements,  we  obtain  computationally  efficient 
and  exact  implementations  of  the  KBf,  reducing  its  computational  complexity  by  a  factor  equal  to  the 
linear  dimension  7  of  the  discretized  grid. 

To  reduce  further  the  computational  complexity,  we  approximate  the  information  matrix  of  the  KBf, 
i.e.,  the  inverse  of  the  error  covariance  matrix,  by  a  sparse  banded  matrix.  Banded  approximations  to 
information  matrices  correspond  to  modeling  the  error  field  in  our  spatial  estimates  at  each  point  in 
time  as  a  reduced  order  Markov  random  field  (Mrf).  It  can  be  shown  that  these  banded  inverse 
approximations  are  optimal  approximations  in  a  maximum  entropy  sense.  We  apply  the  block-banded 
structure  to  derive  several  recursive  algorithms,  which  relate  the  constituent  blocks  in  the  error 
covariance  matrix  to  its  diagonal  and  its  upper  diagonal  entries.  Due  to  the  special  sparse  nature  of  the 
measurements,  most  of  the  off  diagonal  blocks  need  not  be  updated.  We  update  only  the  main  and  the 
upper  diagonal  blocks.  Any  other  block  if  needed  may  be  derived  from  the  Mrf  relationships.  This 
reduces  the  computational  effort  for  the  single  layer  NLOM  by  at  least  0(7  ")  over  the  direct  KBf.  The 
storage  requirements  are  reduced  by  0(7 ).  Our  implementations  provide  additional  efficiencies  with 
multiple  layers. 

WORK  COMPLETED 

Assimilation  of  the  altimetry  data  has  been  investigated  for  the  latest  spherical  version  of  the  NLOM. 
The  derivation  of  the  data  assimilation  algorithm  based  on  the  KBf  was  completed  and  fine-tuned  for 
improved  performance.  We  ran  experiments  with  the  NLOM  at  a  finer  resolution  of  V*  degree  for  a 
more  complex  three  and  a  half  hydrodynamic  NLOM  in  several  regions  of  interest.  Experiments  were 
initially  run  using  both  analytical  winds  and  the  more  realistic  annually  averaged  Rosential  winds. 

In  the  near  future,  extension  to  real  Topex  /  Poseidon  data  is  planned.  This  will  make  our 
experimentation  more  realistic. 

RESULTS 

We  provide  results  for  a  V*  degree  three  and  a  half  layer  hydrodynamic  NLOM  lying  within  -10°  < 
latitude  <  10°,  306°  <  longitude  <  374°.  Figurel  shows  contour  plots  of  the  estimates  of  the  sea  surface 
height  (SSH)  for  the  top  layer  after  5  Topex/  Poseidon  simulated  repeat  cycles.  The  top  plot  in  Figure  1 
is  the  “truth”,  the  middle  plot  in  figure  1  is  the  spin  up  with  no  assimilation,  and  the  bottom  plot  in 
figure  1  is  the  field  after  data  assimilation.  Close  visual  comparisons  of  the  plots  reveal  that  there  are 
several  major  differences  in  the  three  plots.  The  higher  resolution  details  like  eddies  are  completely 
absent  from  the  middle  plot  of  figure  1  while  most  of  them  appear  at  the  correct  position  in  the 
assimilated  plot  at  the  bottom  of  figure  1. 


In  figure  2,  we  plot  for  the  top  two  layers  the  time  evolution  of  the  relative  mean  square  error  (RMSE). 
The  RMSE  is  defined  as  the  ratio  of  the  MSE  of  the  data  assimilated  run  over  the  MSE  of  the 
deterministic  run.  With  time,  the  RMSEs  for  both  layers  drop  suggesting  an  improvement  in  the  data 
assimilated  prediction.  The  improvement  is  higher,  roughly  25%  at  the  end  of  the  experiment,  for  the 
upper  layer  for  which  the  satellite  provides  direct  measurement  data.  No  data  is  assimilated  for  layer  2, 
yet  we  see  an  overall  improvement  of  about  10%.  Our  KBf  implementation  successfully  correlates  the 
two  layers  and  projects  the  surface  information  into  the  lower  layers. 

IMPACT/APPLICATION 

We  have  presented  our  work  at  the  Third  World  Meteorological  Organization  (WMO)  International 
Symposium  on  Assimilation  of  Observations  in  Meteorology  and  Oceanography,  Quebec  City,  Canada, 
see  reference  [1]  below.  We  obtained  good  and  encouraging  reviews  from  the  participants. 

A  full  paper  on  this  work  was  also  accepted  and  will  appear  in  the  November  1999  issue  of  the  IEEE 
Transactions  on  Image  Processing,  see  reference  [2]. 

The  project  is  in  close  collaboration  with  the  research  group  of  Dr.  Hurlburt  at  the  NRL  Stennis  Space 
Center.  Our  goal  is  to  make  available  the  results  to  NRL.  This  work  will  lead  to  improved  ocean 
forecasts  that  can  be  used  with  the  Navy  fleet  nowcast  and  forecast  programs. 

TRANSITIONS 

The  code  will  be  made  available  to  NRL. 
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1:  Contour  plots  of  the  sea  surface  height  (SSH):  Top:  SSHfor  the  actual  truth;  middle:  SSH  with 
no  data  assimilation;  and  bottom:  SSH  with  data  assimilation  from  satellite  altimetry.  Darker 

regions  are  deeper ;  lighter  regions  are  higher. 


2:  RMSE  of  the  SSHfor  the  top  two  layers.  The  lower  curve  is  the  RMSEfor  the  surface  layer 
while  the  upper  curve  is  the  RMSE  for  the  second  lower  layer.  No  data  is  assimilated  for  the 
lower  layers,  only  altimetry  data  for  the  top  layer  is  assimilated  yet  reconstruction  of  all 

layers  is  significantly  improved. 


